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ABSTRACT 

We present a model for the interstellar medium of NGC4214 with the objective to probe the physical conditions in the two main 
star-forming regions and their connection with the star formation activity of the galaxy. We used the spectral synthesis code Cloudy 
to model an Hii region and the associated photodissociation region (PDR) to reproduce the emission of mid- and far-infrared fine- 
structure cooling lines from the Spitzer and Herschel space telescopes for these two regions. Input parameters of the model, such as 
elemental abundances and star formation history, are guided by earlier studies of the galaxy, and we investigated the effect of the mode 
in which star formation takes place (bursty or continuous) on the line emission. Furthermore, we tested the effect of adding pressure 
support with magnetic fields and turbulence on the line predictions. We find that this model can satisfactorily predict (within a factor 
of ~2) all observed lines that originate from the ionized medium ([Siv] lO.Spm, [Neiii] 15.6/jm, [Siii] 18.7yum, [Siii] 33.5yum, and 
[Oiii] 88yum), with the exception of [Neii] 12.8pm and [N ii] 122pm, which may arise from a lower ionization medium. In the PDR, 
the [O i] 63pm, [O i] 145pm, and [C ii] 157pm lines are matched within a factor of ~5 and work better when weak pressure support is 
added to the thermal pressure or when the PDR clouds are placed farther away from the H ii regions and have covering factors lower 
than unity. Our models of the H ii region agree with different evolutionary stages found in previous studies, with a more evolved, 
diffuse central region, and a younger, more compact southern region. However, the local PDR conditions are averaged out on the 
175 pc scales that we probe and do not reflect differences observed in the star formation properties of the two regions. Their increased 
porosity stands out as an intrinsic characteristic of the low-metallicity ISM, with the PDR covering factor tracing the evolution of the 
regions. 

Key words, galaxies: dwarf - galaxies: star formation - galaxies: individual: NGC 4214 - infrared: ISM - techniques: spectroscopic 
- radiative transfer. 


1. Introduction 

The physical state and structure of the interstellar medium (ISM) 
are important parameters for understanding the star formation in 
a galaxy. In a typical star-forming region, young massive stars 
are born and start to illuminate their parental cloud. UV pho¬ 
tons ionize the surrounding medium, creating H ii regions, while 
the transition to the neutral atomic or molecular phase occurs at 
higher visual extinction, where the material is more effectively 
shielded. Far-UV (FUV) photons control the chemical activ- 
ity in these regions, namely the photodissociation region (PDR; 
Tielens & Hollenbachlll985l) . By studying the latter, we can in¬ 
vestigate the conditions of the molecular clouds, which in turn 
will be potential sites for the next episode of star formation. 

How does the propagation of radiation and the ISM compo¬ 
sition affect ISM observables in low-metallicity galaxies? Ad¬ 
dressing this question is important to understand the evolution 
of low-metallicity galaxies, which undergo more bursty star for¬ 
mation than normal galaxies. Nearby star-forming dwarf galax¬ 
ies present distinct observational signatures compared to well- 
studied disk galaxies. Dwarfs are usually metal poor. Hi rich, 
and molecule poor as a result of large-scale photodissociation 
(e.g., iKunth & OstlinI l200(]l : iHunter et alJ 120121 : ISchruba et alJ 


l2012h . Mid-IR (MIR) and far-IR (FIR) observations have re¬ 
vealed bright atomic lines from Hit regions ([ Siii1, [Neml, 


[Neil l, rOiiik etc.) and PDRs (ICiil, [Oij) (e.g.. Hunter et al.l 


20011: iMadden et al.l 120061: IWu et al.l l20o4 iHunt et al.l 120101: 


Cormier et al.l 20151) . Their spectral energy distributions (SEDs) 


are also different from spiral and elliptical galaxies and in¬ 
dicative of altered dust properties, with a relatively low abun¬ 
dance of polycyclic aromatic hydrocarb ons (PAHs) and per¬ 
haps a different dust composition (e.g .. i Maddenet^ 120061: 
iGalliano et al.l 120081: iRemv-Ruver et all 20131) . It is still un¬ 
known, however, whether these differences between dwarf and 
disk galaxies are the direct result of recent star formation ac¬ 
tivity shaping the ISM or instead a consequence of the low- 
metallicity ISM that is independent of star formation activity. 
To answer this, one needs to observe tracers of the interplay 
between the ISM and various stages of star formation activ¬ 
ity. While there are now a number of important studies avail¬ 
able on PDR properties modeling FIR li nes on large scales in 
various extragalactic environments (e.g., 

Wstae^^ 2010l:lGracia-Carpio et al.l201 1 


Kaufman et al. 12006 


Cormier et al.l2012 


20131) or in our Galaxy under solar-metallicity 
.g., Cubick et'^l2008t iBernard-Salas et al1l2012l 


based on observations performed with Spitzer and Herschel. 


Parkin et al. 

condi tions (e.;_ _ 

l2015l) . only a few studies are published on individual extraga lac- 
tic regions (iMookeriea et al.ll2011[ iFebouteiller et al.ll201^ . Of 
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Right ascension 


Fig. 1. Three-color image of NGC4214 using the HST WFC3 
filters F438W (B, blue), F502N ([Out], green), and F657N 
(F1 q'-i-[Nii], red), downloaded from the Hubble Legacy Archive 
(http: //hla. stsci. edu/1. 


particular interest are dwarf galaxies, where the effect due to ra¬ 
diative feedback is expected to be most significant. The goal of 
this paper is to investigate how the low-metallicity ISM reacts 
under the effects of star formation in regions that have undergone 
different histories. The nearby low-metallicity galaxy NGC4214 
provides an excellent environment to perform this experiment 
because it has well-separated star-forming centers, one hosting a 
super star cluster, which allows us to study the effects of extreme 
star-forming conditions on the surrounding ISM. 


N GC4214 is a n e arby irregular galaxy located 3Mpc 


away (iDalcanton et al.l 120091) with a metallicity of ~0.3Zq 
(iKobulnickv & Skillmanlll99CT and a wealth of ancillary data. 
It shows various morphological characteristics su ch as H i holes 
and shells and a spiral pattern (lMcIntvrelll998h . NGC4214 is 
known to host two main, well-defined star-forming regions with 
recent activity (Fig.[T]i. The largest of the two regions is found in 
the center of the galaxy (also referred to as NW or region I) and 
contains several clusters, including a super star cluster, while the 
second region is found to the southeast (also referred to as SE 
or region II) and is younger and more compact. Using near-IR, 
optical, and UV data, several studies have constrained the ages 
of the clusters in the two main regions, which show evidence for 


recen t star formation l Ubeda et al.l 12007 at ISollima et al.l 1201^ 
l2014l) . [^chruba et al.l (l2012h have measured the ongoing SFR of 
NGC4214 to be 0.12Moyr“'. The galaxy seems to have main¬ 
tained its star formation in the past lOGyr at an average rate of 
~0.02Mo yr“\ with a prolonged star formation episode that oc¬ 
curred about 3 Gyr ago and several shorte r bursty events within 
the past Gyr at a rate of 0.05-0.12 Mp yr * (iMcOuinn et al.l2()T()t 
I Williams et al.ll20Tll) . 


In this paper, we present observations of MIR and FIR fine- 
structure cooling lines in NGC 4214, which provide key diagnos¬ 
tics of the physical conditions of the ISM. We focus our analysis 
on the two main star-forming complexes. The line emission is 
analyzed with radiative transfer models to characterize the ISM 
conditions. We take into account directly observed star forma¬ 
tion histories and explore how they affect the IR line emission. 
Photometry is used for the energy budget of the models. The 
structure of this paper is the following: Sect. |2] describes the 


data. Sect. [3] describes the model, and the results are presented 
in Sect. |4] We summarize and discuss our results in Sect. |5] 


2. Data 

2.1. Herschel data 


We used observations of NGC 4214 obtained by the PACS in¬ 
strument (iPoglitsch et alj20lO) onboard the Herschel Space Ob- 
serv atorv (iPilbratt et al. 20101) as part of the Dwarf Galaxy Sur¬ 
vey (iMadden et al. l2013l) . The list of observations can be found 
in Table [U The photometry data at 70/tm, 100/rm, and 160/tm, 
with respective b eam sizes (FWHM) of 5.6 ", 6.7", and 11.3", 
were published bv iRemv-Ruver et alJ (l2013h . These bands cover 
the peak of the SED originating from the reprocessed stellar 
light by the dust. The spectroscopy comprises observations of 
the [Out] 88/tm and [Nii] 122//m lines, which trace the ionized 
gas, as well as the [Cii] 157//m, [Oi] 63/tm, and [Oi] 145/rm 
lines, which trace the PDR. The data consist of small map¬ 
pings of 5 X 5 rasters separated by ~16" for [Oiii] 88/tm and 
[Oi] 63jum and 3x3 rasters separated by ~24" for the other 
lines, ens uring a uniform cover age of 1.6' x 1.6'. Originally pre¬ 
sented in ICormier et akl (120101) . the PACS spectral data were re¬ 
processed wi th the red uction and analysis s oftware HIPE user 
releas e v. 11 (IOttll2010h and PACSman v.3.5 (iLebouteiller et al.l 
l2012h . With the improved calibration and definition of the re- 
gi ons, flux maps are g lobally consistent with those published 
in ICormier et al.l (1201 oh and line ratios agree within 30%. Flux 
maps of the [Cii] 157//m and [Oi] 63/tm lines are shown in 
Fig. 13 The associated error maps include data and line-fitting 
uncertainties, but not calibration uncertainties, which are on the 
order of 15%. The FWHM is 9.5" below 100/tm and 10", 11", 
12" at 122/im, 145/tm, 160jum, respectively. All maps were con¬ 
volved to the [Cii] 157jum resolution of ~12", which at the dis¬ 
tance of NGC 4214 corresponds to a physical scale of 175 pc. In 
both the photometry and spectroscopy dat a sets, the convolutio ns 
were performed using kernels provided bv lAniano et iDdmIE 


2.2. Spitzer data 

NGC 4214 was observe d with the three inst ruments onboard the 
Spitzer space telescope (IWerner et al.ll2004l) . We used the MIPS 
24jum o bservations obtai ned within the Local Vol ume Legacy 
Survey (iDale et al.l l2009l) that were processed by iBendo et al.l 
(l2012l) . The MIPS 24yum map, which has an original FWHM of 
5.9", was convolved to a resolution of ~12" to match that of the 
PACS data. 

The IRS observations (program ID 3177, PI. Skillman) con¬ 
sist of small mapping s of the two main star-forming regions in 
high-resolution mode (iHouck et al.ll2004l) . We extracted the data 
from the Spitzer Heritage Archive (see Table [p and proc essed 
them with the software CUBISM vl.8 (ISmithet al.ll2007h . We 
used the default mapping procedure and bad pixel removal to 
produce spectral cubes with pixel sizes 2.26" for the Short-High 
module and 4.46" for the Long-High module. We then created 
surface brightness maps for all spectral lines of interest - [S iv] 
10.5//m, [Neii] 12.8/tm, [Neiii] 15.6//m, [Siii] 18.7//m, [Siii] 
33.5pn\, which all trace Hii regions - in the following way: for 
each pixel of the cube, we extracted the signal with a range of 
+0.7jum around the line and fit a polynomial of order two for the 
baseline and a Gaussian for the line with the IDF routine mpf it. 
For a more stable fit, the peak of the Gaussian is required to be 

* http://WWW.astro.princeton.edu/-ganiano/Kernels/ 
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Fig. 2. Maps of the [Cii] 157^m, [Oi] 63/jm, and TIR emission in NGC4214. Units are Wsr“'. The two star-forming regions, as defined 
in Sect. 12.41 are outlined with red contours. The right panel shows the Spitzer IRS mapping strategy. Orange: Long-High module coverage; cyan: 
Short-High module coverage; gray background: [Cii] 157yum map. 


Table 1. List of Herschel and Spitzer observations. 


Instrument 

Herschel data 

Observation Identification number (OBSID) 

PACS phot. 

PACS spec. 

1342211803,1342211804,1342211805, 
1342211806 

1342187843,1342187844,1342187845, 
1342188034,1342188035,1342188036 

Instrument 

Spitzer data 

Astronomical Observation Request (AOR) 

MIPS 24pm 

22652672,22652928, 22710528, 22710784 


IRS Short-High 

On-source: 10426368,10426624, 10426880, 10427136 

Background: 13728256,13729792, 13730304, 13730816, 
13733120,13762304,13767424,13767936, 
13768448,13768704,13769728,13770496, 
13773568 

IRS Long-High 

On-source: 10424832,10425088, 10425344, 10425600, 

10425856,10426112, 10427392 

Background: 13728768,13729792, 13730304, 13730816, 
13733120,13763328,13764352,13765376, 
13767424,13767936,13768448,13768704, 
13769728,13770496,13773568 


positive, the position of the peak is expected within one instru¬ 
mental FWHM of the rest wavelength, and the width is limited 
to the instrumental resolution {R - 600). Finally, we added ran¬ 
dom noise and iterated the fit 300 times to estimate the best-fit 
parameters as the median of the resulting parameters and the er¬ 
ror on those parameters as the standard deviation. Error maps 
again include data and line-fitting uncertainties, but not calibra¬ 
tion uncertainties, which are on the order of 5%. The coverage 
of the star-forming regions of the galaxy in the IRS maps is only 
partial, as shown in Fig.|2] No integrated values for the flux of 
the whole regions could be retrieved. To obtain a representative 
value for the line flux in each region, we regridded the IRS maps 
to that of the [Oiii] 88//m map. Then we selected the pixels that 
appear in both maps and scaled the emission of these pixels to 
the [O III] 88;um line to infer the corresponding line fluxes for the 
star-forming regions as a whole. 


We focus on these selected IRS lines because they are among 
the brightest MIR fine-structure cooling lines and can be used 
as reliable diagnostics of the physical conditions in Hii re¬ 
gions. In general, the intensity or luminosity ratio of two lines 
of the same element but different ionization level is indicative 
of the radiation held hardness. Such diagnostics are the [Neiii] 
15.6iU m/rNeiil 12.8iUmo rthe [Siv] 10.5//m/[S m] 18.7/im ratios 
(e.g.. IVerma et al.l l2003h . which are insensitive to the density 
because of their high critical densities (see Table la. Accord¬ 
ingly, species of the same ionization level but different transition 
are indicative of the electron density as a result of the differ¬ 
ent cr itical densities for each transition (lOsterbrock & Ferlandl 
l2006h . Examples are the [Siii] 18.7//m/[S m] 2>3.5pm, [Neiii] 
15.6iU m/rNeiiil 36.0Atm , or [Nii] 122jum/[Nii] 205jum ratios 
(e.g., iRubin et al.l H994ll . These diagnostics are insensitive to 
the temperature inside the Hii region. Unfortunately, the [Neiii] 
36.0//m and [N ii] 205//m lines fall at the edge of the IRS and 
PACS wavelength ranges, respectively, where the spectra are too 
noisy to detect or derive a reliable line ratio for the two star¬ 
forming regions. Therefore we relied on the [S iii] line to probe 
the electron density. 

2.3. Total infrared luminosity map 

To construct a total infrared (TIR) luminosity map of the galaxy, 
we combined t he MIPS 24juni and PACS 70, 100, and 160jum 
data, following iGalametz et al.l (l2013h : 



We used the values of the coefficients, c;, from Table 3 of their 
paper: [c 24 ,C 70 ,cioo,cieol = [2.064,0.539,0.277,0.938]. This 
method, although slightly less accurate than a direct integration 
of a well-sampled SED, does not require degrading the resolu¬ 
tion of our data beyond the PACS 160//m beam and is sufficient 
for our modeling purposes to estimate the energy budget in the 
star-forming regions. The Ltir map is shown in Fig.|2 

2.4. Defining the star-forming regions 

To define the apertures for the main star-forming regions, we 
set a threshold for the signal-to-noise ratio (S/N) equal to 5 in 
each individual PACS IQpm, lOOjum, and 160/tm photometry 
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and PACS spectral map. We masked all pixels below this S/N and 
drew the contours, which include all the remaining unmasked 
pixels separately for the photometry and the spectroscopy maps. 
Because the emission in the photometry maps is more extended, 
we kept the contours from the photometry and used these aper¬ 
tures throughout the analysis to define the two star-forming re¬ 
gions, as shown in Fig.|2] This means that pixels in the spec¬ 
troscopy maps that are below the S/N threshold but within the 
region contours are still counted. The fluxes and uncertainties 
for the line and TIR emission were measured taking into account 
all pixels in each region. They are reported in Table |2] The ISM 
emission (gas and dust) peaks in these two regions, and most of 
the line fluxes are twice as high in the central region as in the 
southern region, except for [Ne iii] 15.6;um (factor 1.3) and [S iv] 
10.5//m, which have lower fluxes toward the central region. This 
hints at different physical conditions in the two regions, which 
we investigate with radiative transfer models. 

3. Description of the modei 

3.1. Model geometry and strategy 

Our objective is to characterize the physical conditions of the 
ISM phases from which the IR emission arises in NGC4214. 
To that end, we u sed the spectral synt hesis code Cloudy v.13, 
last described by iFerland et al.l (l2013h . We performed a mul¬ 
tiphase detailed modeling of the ISM for which we combined 
line and continuum e mission, following the method described in 
ICormier et al.l (l2012h . We considered the two main star-forming 
regions of NGC4214: the most evolved central region (NW-I) 
and the southern region (SE-Il). Here we present the main as¬ 
pects of the model and how it is applied to each region. The 
components/lSM phases of the model are 

1. a central source of radiation, 

2. an ionized medium component (Hit region) that surrounds 

the central source, 

3. a neutral medium (PDR) surrounding the H ii region. 

This method assumes a single radiation source responsible for 
the observed SED of the studied region. In other words, we 
took all of the different sources (star clusters) and the surround¬ 
ing clouds from which they have formed and represented them 
with one central source and one surrounding cloud. We thus tar¬ 
geted the integrated properties of each region. In practice, we 
mixed components of the ISM that have different composition 
and properties and blend them in a single system. The applied 
geometry is spherical. The source is in the center and the illumi¬ 
nated face of the cloud lies at a certain distance that we call inner 
radius. In our case, the effective geometry is one-dimensional 
plane-parallel because the cloud forms a thin shell and its dis¬ 
tance from the radiation source is large. 

The radiation source, representative of the stars that popu¬ 
late the clusters of the star-forming region, illuminates a cloud 
of dust and gas. It controls the ionization parameter, U, which 
characterizes the held and is defined as the ratio of the incident 
ionizing photon density to the hydrogen density. Hard UV pho¬ 
tons from the source ionize hydrogen and form the H ii region. 
As this radiation is transmitted through the cloud, it is attenu¬ 
ated and thus becomes softer, which decreases its influence on 
ionization. However, it still controls the processes further in the 
cloud (in the PDR). 

The adopted strategy is to treat the H ii region first and then 
use the H ii region parameters as input for t he PDR modeling . 
This allows for a self-consistent approach (lAbel et al.l l2005h . 


whi ch is usually not d irectly available in standard PDR codes 
(see iRollig et^l2007l for a comparison of PDR codes), and is 
important to accurately derive the radiation held that impinges 
on the PDR. We first ran the simulation until the end of the H ii 
region, choosing to stop the simulation when the (electron) tem¬ 
perature reached 500 K to ensure that the model had transitioned 
to the atomic phase. We iterated to optimize our parameters so 
that they matched the observed emission of known H ii region- 
diagnostic lines ([Siv] 10.5/rm, [Neii] 12.8/rm, [Neiii] 15.6/rm, 
[Sill] 18.7/rm, [Siii] 33.5/rm, [Out] 88/rm, and [Nii] 122/um). 
Then we fed the result of this model to the PDR and compared 
the predictions to the remaining three PDR lines observed: [O i] 
63//m, [Oi] 145/tm, and [Cii] 157/rm, choosing a visual extinc¬ 
tion of 10 mag as the stopping criterion. At this point, the gas 
temperature had fallen to roughly 10 K. 

3.2. Model parameters 

We constrained the properties of the star-forming regions by 
varying some of the parameters that control the physics of the 
models while keeping others fixed. The main parameters that we 
consider are 

1. a radiation field source: shape, age, luminosity (varied), 

2. the hydrogen density of the ISM, tin (varied), 

3. the ISM gas elemental abundances (fixed), 

4. the inner radius, rin (varied), 

5. the magnetic held, B (fixed), 

6. the turbulent velocity, Uturb (fixed). 

Parameters that are fixed were set to values from the literature. 
The other parameters were varied inside a range whose width 
reflects the dispersion of published measurements or of the data. 
The main parameter of interest for this study is the radiation 
held, which was varied within a range guided by studies of the 
star formation history. 

3.2.1. Hydrogen density («h) 

We performed our simulations assuming pressure equilibrium. 
As the model proceeds through consecutive zones of the cloud, 
it keeps the pressure constant. Thus, the density of the medium 
varies to satisfy this equilibrium. The initial density that we 
specified in the models is the density at the illuminated face 
of the cloud, where the Hit region starts. The initial values 
and the range we probed are motivated by the observed [S iii] 
18.7pm/[Sra] 33.5/rm ratio in the Hit region, which is known 
to be sensitive to t he electron density in the range 10^-10"^ cm“^ 
(iRubin et al.ll994h . We therefore let the initial density vary in the 
range 100-300 cm“^ for the central region and 300-600 cm“^ for 
the southern region with a common step of 25 cm“^. By iterating 
this procedure, we constrained the best values for the density at 
the beginning of the H ii region. 

3.2.2. Inner radius (rin) 

In our spherical geometry, the source is at the center and is sur¬ 
rounded by a cloud. The illuminated face of the cloud lies at a 
certain distance rin. This is not a strictly physically constrained 
parameter because the setup we used does not realistically model 
each cluster, but instead tries to mimic a whole region and repro¬ 
duce its emission. The variation of this radius changes the pho¬ 
tons flux and thus is expected to affect our results. We let ri„ vary 
from 1 to 100 pc for both regions. 
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Table 2. Observed MIR and FIR fluxes for the line and broadband emission. Uncertainties on the fluxes include data and line-fitting uncertainties, 
but not calibration uncertainties, which are on the orde r of 5% for the Spitzer lines and 15% for the Herschel lines. Critical density and ionization 
potential values are taken from ICormier et aP l l2012f) . Critical densities are noted [e] for collisions with electrons and [H] for collisions with 
hydrogen atoms. 


Line 

Wavelength 

(jim) 

Elux + uncertainty 
(xlO^*® Wm-2) 

Region I Region II 

Critical density 
(cm“^) 

Ionization 

potential 

(eV) 

[Siv] 

10.51 

5.68 + 0.21 

8.40 + 0.10 

5 X lO'* [e] 

34.79 

[Ne II] 

12.81 

8.98 + 0.22 

4.13+0.11 

7 X 10^ [e] 

21.56 

[Ne III] 

15.56 

18.70 + 0.14 

14.25 + 0.08 

3 X 10^ [e] 

40.96 

[Sm] 

18.71 

11.78 + 0.20 

6.85 + 0.08 

2 X 10“ [e] 

23.34 

[Sm] 

33.48 

18.71 + 0.27 

8.20 + 0.12 

7 X 10^ [e] 

23.34 

[Oi] 

63.18 

10.11+0.35 

4.06 + 0.21 

5 X 10^ [H] 

- 

[Oiii] 

88.36 

31.86 + 0.62 

13.50 + 0.40 

5 X 10^ [e] 

35.12 

[Nil] 

121.90 

0.44 + 0.20 

0.14 + 0.08 

3 X 10^ [e] 

14.53 

[Oi] 

145.52 

0.65 + 0.09 

0.32 + 0.07 

1 X 10^ [H] 

- 

[Cii] 

157.74 

26.34 + 0.33 

10.05 + 0.21 

50 [e], 3 X 10^ [H] 

11.26 


Broadband 

Wavelength 

(pm) 

Elux density + uncertainty 

(Jy) 

Region I Region II 

MIPS 

24 

0.67 + 0.01 

0.48 + 0.01 

PACS 

70 

7.36 + 0.21 

3.72 + 0.12 

PACS 

100 

7.91+0.19 

4.10 + 0.11 

PACS 

160 

6.07 + 0.10 

3.23 + 0.06 

Ltir (ergs^^) 

3-1100 

5.25 X 10^^* 

3.02 X lO'^* 


Table 3. Elemental abundances in NCC4214. Values (in logarith mic 
scale) for NCC 4214 are taken fr omiKobulnickv & SkillrnarJ d 1 996h and 
solar values from lAsnlund et akl d2009h . 


Abundance 

Region I 

Region II 

Solar value 

[O/H] 

-3.795 + 0.05 

-3.64 + 0.04 

-3.31+0.05 

[S/H] 

-5.380 + 0.06 

-5.21+0.06 

-4.88 + 0.03 

[N/H] 

-5.094+ 1.00 

-5.02 + 0.10 

-4.17 + 0.05 

[Ne/H] 

-4.535 + 0.11 

-4.51+0.08 

-4.07 + 0.10 

[C/H] 

-4.295 + 0.30 

-4.14 + 0.30 

-3.57 + 0.05 


3.2.3. Elemental abundances 


star fo rmat ion history of the gala xy presented in lMcOuinn et al.l 
(1201 nh and lWilliams (1201 111 . we tested the following cases: 

- For the central region, we considered two limiting scenar¬ 
ios: (i) a single-burst star formation event and (ii) a con¬ 
tinuous star formation model, with SFR=O.O7M0 yr^'. The 
ages of the clusters were varied within a range of (i) 1- 
20 Myr (in steps of 0.5 Myr) and (ii) 200-1000 Myr (in steps 
of200Myr). 

- For the southern region, we considered a single-burst event 
with an age that varied from 1 to 20 Myr (in steps of 
0.5 Myr). 


Elemental abundances in the models were set to the observed 
values for oxygen, su l fur, ni trogen, and neon, taken from 
iKobulnickv & SkillmanI (Il996l) . Some measurements partially 
cover our defined regions, and we adopted them as representa¬ 
tive. Exclusively for carbon, we scaled its abundance according 
to the study on the de pendence of log(C/(9) on metallicity by 
lizotov & ThuanI (Il999l) . Eor other elements, we used the default 
ISM composition of Cloudy and scaled the abundances to our 
metallicity (1/3). The values used are indicated in Table[3] 


3.2.4. Radiation field - shape and energy 

We used the code Starburst99 (iLeitherer et alJ20T(ih to produce 
a stellar population spectrum that serves as input for our mod¬ 
els. We chose a Kroup a initial mass function between 0.0 8 and 
120Mo (iKrounalbOOlh . as done in [Andrews et al.l (12013h . and 
Padova asymptotic giant branch tracks with Z=0.004. As dis¬ 
cussed above, we did not model each cluster individually, but we 
used integrated emission from the entire star-forming regions in¬ 
stead. We tried to be as close to the shape and intensity of the 
radiation field of the Hii regions as possible. Motivated by the 


A fix ed mass of 10^ Mq (typical cluster mass in ISollima et al.l 
1201 3h was considered for the single bursts, where the stars 
are created at once (delta burst). Our value for the SER 
(0.07Moyr^') is representative of the ‘average’ rate at which 
this galaxy formed stars within the past 1 Gyr of its history. The 
ages of the clusters in both regions were guided by values from 
[Ubedaet^alJj 2007a), Sollima et al.l (1201 3h . and Sollimaet al.l 
( 2014 ). Ubeda et al. (l2007a ) found 2-7Myr for region I, along 
with extended clusters in the same region with ages of 150 to 
190 Myr. In r e gion II, they found ages spreading around 2 Myr. 
ISollima et akl (120131 1201 4l) reported a larger spread in ages. In 
region I, they obtained ages around a median of 14 Myr, and for 
the more extended clusters the ages lie between 10 and 300 Myr. 
Eor region II, they found a median age of ~20 Myr. 

Eor the luminosity emitted from each region, we chose to use 
the TIR luminosity as a first approximation of the luminosity of 
the starburst. This choice implies the assumption that all radia¬ 
tion from the clusters in the region is reprocessed by the dust and 
thus emitted at longer wavelengths. In doing so, we kept in mind 
that there can be processes that we did not model (UV escape 
fraction or a diffuse ionized medium, for example) and that can 
contribute to this radiated energy (see Sect.|3. 
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3.2.5. Magnetic field strength (B) and turbulent velocity (wturb) 

Magnetic fields and t urbulence play an importa nt role in struc¬ 
turing the ISM (e.g.. iMcKee & Ostrikedl^OTh . When Cloudy 
solves the pressure equilibrium for each zone of the modeled 

cloud, a magnetic pressure term equal to Pb = — is included in 

Stt 

the equation of state along with a turbulent pressure term equal 

to = 2.8 ■ 10® ■ 3 ■ { i )" [cm-^K], for 

10^ cm^^ lkms“‘ 

isotropic turbulent motions, where «h is the total hydrogen den¬ 
sity and Cturb is the turbulent velocity (see the Hazy documenta¬ 
tion of Cloudy for more information). 

The magnetic field of NGC4214 was measured by 
iKenlev et al.l (1201 ih using multiwavelength radio emission. The 
reported field strength in the center of the galaxy is 30 pG, and 
the pressure term due to this field has the same order of magni¬ 
tude as the hot gas and the gravitational contributions. Since it is 
not well known how the observed magnetic field might affect our 
observed line intensities, we excluded it from our default models 
and tested one case with a magnetic field strength of 30pG. 

Another potential energy source to consider can arise from 
the dissipation of turbulence. Turbulent energy is converted into 
thermal energy as it cascades from large scales to small scales 
through dissipation. However, we did not resolve size scales 
for which we can measure this. Cloudy does not attempt to 
model the dissipation mechanism, but assumes a simple ther¬ 
mal energy source based on line width. The turbulent velocity 
was set to a value of 1.5kms“^ by default in our models, and 
we tested two other cases: one case with an intermediate tur¬ 
bulent velocity (iiturb=3kms“* or FWHM=5kms“*) that corre¬ 
sponds to the approximate line width observed in the CO(l-O) 
data by IWalter et al.l (1200 ih . and one case with a high turbulent 
velocity (t;iiirb-50kms~^) as found in the diffuse ionized gas b y 
IWilcots & Thurovil (1200 ih and used also in lKenlev et al.l (1201 ih . 

Nevertheless, we explore the effects of excluding or includ¬ 
ing magnetic fields and turbulence in Sect. 14.21 

3.3. Determination of best-fitting modeis 

We aim to converge on a unique parameter set that best describes 
the conditions of the regions. We first ran models for which we 
varied the parameters in a coarse grid to narrow down the param¬ 
eter space, using ranges of values found in the literature to start 
with. We then used the optimization option of Cloudy, which au¬ 
tomatically varies the specified parameters in a finer grid to find 
the optimal solution. 

We computed the average denoted for each model by 
comparing the observed fluxes of [S iv] 10.5;um, [Ne iii] 15.6/tm, 
[S III] 18.7/im, [S III] 33.5pm, and [Oiii] 88//m to the fluxes pre¬ 
dicted from the radiative transfer calculation. These lines are the 
most luminous and most strongly correlated with the Hii re¬ 
gion (as opposed to [N ii] and [Ne ii] which, from experience, 
can arise from other phases). We refer to these as the optimized 
lines in Table|5] The goodness of the line emission fit is given by 
low^^ values. The optimization method of Cloudy searches the 
minimum of^^ that is defined as 

v2=lyy2 = ly _ ~ _, (2) 

n 4—1 ' n 4lj (jnin{0p, Mi] X cr;)^’ 

where n is the number of lines optimized and x] are the y val¬ 
ues of the individual optimized lines. M, and O, are the mod¬ 
eled and observed fluxes, and cr, is the fractional error on the 
observed flux (uncertainty/flux) with calibration uncertainties 


added in quadrature to the measured uncertainties described in 
Sect. 12 We have five observables (the ionic lines listed above) 
and varied four parameters (cluster age, source luminosity, hy¬ 
drogen density, and inner radius). 


4. Results 

In this section we present the model results for the two regions 
according to the star formation histories considered. Parameters 
of the best-fitting models and their corresponding values are 
reported in Tables|4]and|5] respectively. 


4. 1. Line emission 


4.1.1. Central region (I): The single-burst mocdel 

For the single-burst star formation event, the best-fitting model 
of the central region has the following parameters: burst age of 
4.1 Myr with luminosity 1.1 x 10"^^ergs"', density of 155cm“^, 
and Lin - 85 pc. The corresponding ionization parameter at the 
illuminating face of the cloud is log(t/) = -2.7. For the Hii re¬ 
gion, all optimized lines are matched within +30%, and [Neii] 
12.8^(m and [Nii] 122yum are underpredicted by a factor of ~8 
and 3, respectively (see Fig. [3] blue bars in the top panel). In 
the PDR, the [Oi] \A5pm and 63pm lines are overpredicted by 
a factor of ~2.5, and the [Cii] 157pm line is matched within 
20%. We note that [Cii] emission can arise from both the neu¬ 
tral and the ionized phases of the ISM, with a potentially non- 
negligible co ntribution fro m the warm ionized medium in the 
Milky Way (iHeilesI [To^ . The Hii region of our best-fitting 
model contributes negligibly to the predicted [Cii] 157pm and 
[Oi] emission (<3%). The ultraviolet radiation fi eld strength Gn 
inside the PDR is 455 in units of the equivalent iHabina (1 1968!) 
flux (1 Go = 1.6 X 10“^ ergcm^^ s '). The input luminosity re¬ 
quired by the model, which all comes out as Ltir in the PDR, is 
about twice as high as the observed Ljir in that region. 

These results represent the simplifying case that the recent 
starburst dominates the star formation in this central region, so 
the line emission can be explained with a single burst. The age 
of the burst in the model nicely agrees with the range of ages 
from Ubeda et al. ( 2007^ and is at the younger end of ages from 
LSollima et al.1 (1201 31 12014 . 


4.1.2. Central region (I): The continuous star formation model 

In the continuous star formation scenario, the best-fitting model 
of the central region has the following parameters: stellar age of 
440Myr with luminosity 1.9 x 10"'^ erg s“', density of 180cm“^, 
and Lin - 62pc. The ionization parameter is \og{U) - -2.5. 
Line predictions for the H ii region and the PDR from this model 
solution are shown in Fig. [3 (gray bars). For all optimized lines 
([Oiii] 88pm, [Siii] 18.7 and 33.5pm, [Siv] 10.5pm, and [Neiii] 
15.6pm), the model matches the observations within ±20%. The 
two other ionic lines, [Ne ii] 12.8pm and [N ii] 122pm, are under¬ 
predicted by a factor of ~1 and 4. In the PDR, [C ii] and the [O i] 
lines are matched within a factor of ~2. The contribution of the 
Hii region to the predicted PDR emission is only 1%. We find 
Go - 1.2 X 10^, which is higher than in the single-burst case. The 
luminosity of the model exceeds the observed Ltir by a factor of 
3.5. 

The results represent the simplifying case of continuous star 
formation dominating this region, with the starbursts being em¬ 
bedded in it. The age of the model agrees well with the star 
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Table 4. Input parameters for the best-fitting models. Values are given at the illuminated face of the cloud. Values in parenthesis for the magnetic 
field strength and turbulent velocities are tested in Sect. 14.21 


Parameter 

Central region 

Single burst Continuous 

Southern region 
Single burst 

Density «h [cm“^] 

155 

180 

440 

Inner radius Hn [pc] 

85.4 

62.1 

22.3 

Stellar age t [Myr] 

4.1 

440 

3.9 

Total luminosity L [erg s“*] 

1.15 X 10^2 

1.91 X lO'^^ 

5.25 X 10^' 

Magnetic field B [pG] 

- 

-(30) 

- 

Turbulent velocity Uturb [km s“* ] 

1.5 

1.5 (3,50) 

1.5 


40 

§ 39 
.03. 

_] 

c5 38 

37 


South region 


39 

CD 

“ 38 

o 

O 

37 


36 

Fig. 3. Results for the central region (top panel) and for the southern region (bottom panel): line emission for the Hit region (left side) and the 
PDR (right side), assuming pressure equilibrium. Green bars represent the observations, blue bars our single-burst model predictions, and gray 
bars with a dashed outline the continuous star formation model predictions. 



Central region 



formation even t in th e window 400-500 Myr ago reported by 
iMcOuinn et al.l (l2010l) . 


4.1.3. Southern region (II) 

The best-fitting model for the southern region is characterized 
by a burst age of 3.9 Myr with luminosity 5.3 x 10"^' erg s“', den¬ 
sity of 440cm~^, and rin - 22pc. The ionization parameter is 
log(t/) = -2.3. Line predictions for the H ii region and the PDR 
are shown in the bottom panel of Fig. [2 This burst found for the 
southern region is slightly younger than the burst in the central 
region, in agreement with llJbeda et alJ (12007 al) . In the Hit re¬ 
gion, the [Out] 88//m, [S iii] 18.7 and 33.5jum, and [S iv] 10.5jum 
lines are reproduced within 30%, while [Neiii] 15.6;um, [Neii] 
12.8//m, and [N ii] 122/im are underpredicted by a factor of 1.7, 
10, and 6, respectively. Feeding this model to the PDR, the [C ii] 
157yum line is underpredicted by a factor of 3.4, while the [Oi] 
lines are both overpredicted by a factor of ~2. The contribution 
of the Hit region to the PDR line emission is only 1%. Go is 
found to be about 3.2 x 10^, which is higher than in the central 


region. The luminosity of the model is 1.7 times higher than the 
observed Ltir in this region. 

4.1.4. Comparison to empirical line ratios 

Physical conditions in the H ii region are mainly determined by 
tracers of the radiation field strength ([Ne iii]/[Ne ii], [S iv]/[S iii] 
18.7//m) and of density ([S iii] 18.7//m/[S iii] 33.5pm). We com¬ 
pare these well-known diagnostic ratios in the two star-forming 
regions in Table|6] In the southern region, ratios of [Ne m]/[Ne ii] 
and [S iv]/[S m] 18.7pm are observed to be about twice as high 
and [Sill] 18.7pm/[S III] 33.5pm marginally higher than in the 
central region, indicating that the radiation field is harder and the 
medium denser. This is indeed what we recover with our best¬ 
fitting models (Table HI, as they match the sulfur ratios well. 

4.2. Effects of magnetic fields and turbulence on the PDR 
temperature and density 

After determining the best parameters for the H ii regions, we ex¬ 
plored the effect of cloud density on the PDR emission in more 
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Table 6. Observed and predicted MIR line ratios for the two star-forming regions. The [Ne iii]/[Ne ii] and [Siv]/[Siii] line ratios are indicative of 
the radiation field strength, and the ratio of the two [S iii] lines is a density diagnostic. 

Ratio Central region Southern region 



Observed 

Burst 

Continuous 

Observed 

Burst 

[Ne iii]/[Ne ii] 

2.082 

13.72 

12.34 

3.448 

28.955 

[Siv]/[S III] 18.7 

0.481 

0.479 

0.477 

1.227 

1.044 

[Siii]18.7/[Siii]33 

0.630 

0.593 

0.601 

0.835 

0.846 


Table 5.x^ values for the two star-forming regions. 

HII region 

Central region 

Southern region 


Burst 

Continuous 

Burst 

Individualvalues: 




[Oiii] 88/im 

0.01 

0.03 

0.01 

[Nil] 122yum 

12.58 

44.63 

82.56 

[Sill] 18.7yum 

2.24 

2.38 

37.69 

[S III] 33.5yum 

8.10 

6.99 

31.89 

[Siv] 10.5|im 

0.76 

0.70 

4.94 

[Neii] 12.8//m 

14442.2 

13842.6 

58101.7 

[Neiii] 15.6/rm 

8.14 

24.65 

224.51 

X~ (all ionic lines) 

2067.72 

1988.85 

8354.76 

(optimized lines) 

3.85 

6.95 

59.81 

PDR 

Central region 

Southern region 


Burst 

Continuous 

Burst 

Individualvalues: 




[Cii] 157//m 

1.57 

31.83 

258.63 

[O i] 63pm 

149.78 

57.05 

85.37 

[Oi] 145/im 

34.14 

11.83 

21.27 

(all PDR lines) 

61.83 

33.57 

121.76 


detail. Density is of critical importance in the emission output 
of the simulation because of the different critical densities of the 
observed lines. Density values quoted so far are representative 
of the Hit region and evolve inside the modeled cloud. Figure |4] 
shows hydrogen density profiles in the clouds for each case pre¬ 
sented in Sect. 14.11 The density starts at the initial value we set 
for each model in the Hit region, remaining practically at the 
same level throughout the H ii region. At the interface between 
the HII region and the PDR, there is a jump in density required to 
keep the model in pressure equilibrium. When pressure is only 
determined by the gas pressure (Pgas = i.e., magnetic 

fields and turbulence are omitted), the temperature difference at 
the phase transition is balanced by a rise in density of 2-3 orders 
of magnitude. Within this frame, the effects of magnetic fields or 
turbulence, implemented as pressure terms in Cloudy, can be un¬ 
derstood as follows: when total pressure equilibrium is assumed, 
they give more support at the phase transition, thus preventing a 
large difference in density between the H ii region and the PDR 
and moderating the density increase at large optical depths. For 
the model presented in Sect. 14.11 where a low turbulence value 
is included but no magnetic fields, representative PDR densities 
are 2 X lO'^cm^^ in the central region for both the single-burst 
and continuous star formation models, and 7 x 10"^cm“^ in the 
southern region (see Fig. Ell. 

The effects of magnetic field and turbulence on the cloud 
density, temperature, and line emission are shown in Fig.|5] We 
present a set of runs for our best-fitting model in the central re¬ 
gion single-burst case (note that we recover similar behaviors for 
the central continuous and southern single-burst cases, as shown 
in Appendix lAli with the magnetic field (B=30//G) and/or the 
turbulence pressure (wturb^l-S, 3, 50kms“') terms switched on. 



0.001 0.010 0.100 1.000 10.0 
Av [mag] 


Fig. 4. Density profiles in the modeled clouds for the central and south¬ 
ern region, which include a turbulence pressure term (i;turb=1 -5 km s“‘). 
Note that the x-axis is logarithmic, so the H ii region (with a constant 
low density) occupies a thin layer of the cloud, stopping at low visual 
extinction (~0.1 mag). 

These terms have no impact on the ionic line emission because 
thermal pressure dominates the pressure balance in the H ii re¬ 
gion, but they noticeably change the emission of the PDR lines. 
When only thermal pressure is considered, the gas density jumps 
to values > 3 x lO'^cm^^ in the PDR and the [Oi] lines are 
overpredicted by an order of magnitude (black bars and dot¬ 
ted line). With only the magnetic field on, all three PDR lines 
are underpredicted by a factor of ~3 (orange bars), but their ra¬ 
tios are kept in the range observed thanks to the lower densities 
achieved (2 x lO^cm^^^). Comparing models with different tur¬ 
bulent velocities, we see that the [C ii] line is best matched for 
low or intermediate velocities because of their moderate den¬ 
sities (~ lO^cm^^) and slightly lower PDR temperatures (at 
Ay ^1-3 mag). Increasing the turbulent velocity reduces the 
predicted [Oi] emission and PDR density. The high-turbulence 
model performs poorly because it has the most dramatic effect 
on the density and line emission. 

To summarize, we find that the best case that simultaneously 
matches all three PDR lines in the central region is the model 
with intermediate turbulent velocity (i;turb=3kms“'), which has 
a density of 8 x 10^ cm“^, but we stress that the main effect of 
turbulence is to reduce the PDR density. 

4.3. Input spectra and SED 

Figure |6] shows the input and output SEDs of the models for the 
two regions. The input SED is the stellar spectrum of the illumi¬ 
nating source modeled with Starburst99 and also includes the 
CMB at millimeter wavelengths. In the central region, the stellar 
spectrum has a wider distribution for the continuous star forma¬ 
tion model (red curve) than the single-burst star formation model 
(black curve), and it is more luminous in the near-IR regime due 
to the presence of old stars. Compared to observations, all input 
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Fig. 5. Effect of turbulence and magnetic fields on the predicted line intensities (top panel), density, and temperature (bottom panels) in the modeled 
cloud for the central region single-burst case. Green bars: observations. Gray bars and solid lines: only low turbulence switched on (t)tu,b=1.5 kms^’, 
default model). Black bars and dotted lines: no magnetic fields and no turbulence. Orange bars and dashed lines: only magnetic field switched on 
(B=30yuG). Red bars and dash-dotted line: only moderate turbulence switched on (!;turb=3 km s“*). Cyan bars and triple-dotted-dashed line: only 
high turbulence switched on (Dturb=50km s“*). 


SEDs fall above the GALEX EUV data because they are unat¬ 
tenuated, and the single-burst input SEDs fall below the 2MASS 
data because they lack old stars. In the central region, the input 
SED of the continuous model, on the other hand, agrees well 
with the 2MASS data. Eor better agreement with the EUV data, 
we estimated the level of extinction required to attenuate the in¬ 
put SEDs. We considered average extinction values E{B - V) of 
0.1, 0.05, and 0.1 mag for the central continuous central single¬ 
burst and southern single-burst models, respectively (dotted lines 
in Eig.lp , which are in the range of values found bv llJbeda et alJ 

(I2007hr 

Eocusing on the output SEDs of the models for the central re¬ 
gion, it is anticipated that the level of the EIR continuum is differ¬ 
ent. The higher the input luminosity of the source, the higher the 
peak of the output SED. Moreover, as the dust temperature rises, 
the peak of the output SED is expected to shift to shorter wave¬ 
lengths. Eor the continuous model, the higher EUV luminosity 
therefore provides more dust heating, explaining the slight shift 
of its peak to shorter wavelengths. We can compare the output 
SEDs to the observed PACS lO/im, lOO/rm, and 160//m fluxes. 
The models agree relatively well with observations for the south¬ 
ern region, but overpredict the EIR continuum emission in the 
central region. This is not surprising because we used as input 
luminosity a higher value than the observed TIR flux and the 
modeled PDR has a high Ay. Better agreement with continuum 
observations requires a model that predicts a TIR flux lower by 


a factor of 2-3, for example, by reducing the covering factor of 
the PDR. We return to this in Sect. 15.1.3] 

5. Discussion 

We have presented models for the two star-forming regions of 
NGC 4214 that work for most of the observed MIR and EIR 
lines. Some discrepancies remain between our models and ob¬ 
servations ([Ne It], [N It], [C ii]/[0 1 ], and Ltir); these are not due 
to the choice of parameter space but rather to missing physics or 
components in our models. In this section, we discuss various as¬ 
pects of our analysis; 1) the discrepancies with observations, and 
we give clues to improve our models, 2) which star formation 
scenario describes the data better, and 3) how our ISM results 
relate to the known properties and evolution of the two regions. 

5.1. Discrepancies between models and observations 

5.1.1. Ionic lines 

In our best-fitting models of the Hit region, the [Nen] 12.8//m 
and [N ii] 122yum lines are systematically underpredicted. These 
lines have lower excitation potentials (21.56 and 14.53 eV, re¬ 
spectively) than the other, better matched ionic lines, and can 
therefore be partially excited outside of the main H ii region. This 
effect can be significant in NGC 4214 because of the poor spatial 
resolution. 
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Fig. 6. Spectral energy distributions of the two star-forming regions: central region single-burst case (left panel), central region continuous case 
(middle panel), and southern region single-burst model (right panel). The black and blue curves correspond to the input and output SEDs, respec¬ 
tively. The dotted lines are the attenuated input SEDs. The data points are the photometry measurements from GALEX EUV, 2MASS J, H, K 
bands, MIPS 24/rm, and PACS at 70yum, lOOyum, and 160pm. In panels for the central region, the dashed curves are scaled versions of the output 
SEDs, considering a covering factor of 0.5 for the PDR. 


Discrepancies between the observed and modeled 
[Ne iii]/[Ne ii] ratio, with the same amplitude as we found 
for NGC 4214, were reported bv iMartm-Hernandez et al.l (l2002h 
for Hit regions observed by ISO/SWS (see their figure 2). 
These could again be related to a mixture of physical conditions 
within the ISO beam, which is also relatively large. For the 
starburst galaxy Haro 11, we have examined the effect of an 
additional low-ionization component (star of effective tem¬ 
perature 35 000 K and hh - 10'“^cm~^). This reproduced the 
observed [Neii] 12.8/um and [Nii] 122 pm emission withou t 
significantly affecting the other ionic lines (ICormier et al.ll20f^ . 
Alternatively, given that the neon lines have the largest energy 
difference, they are more sensitive to the underlying stellar 
atmosphere models than the sulfur lines. Constraining those 
models is beyond the scope of this paper, and we relied on the 
sulfur lines as being more robust diagnostics of the H ii region 
conditions in NGC4214 (i.e., models presented in Sect. 14.lb . 

We also assessed the effect of including the [Nett] 12.8/t(m 
line in the best-fitting solution-tracking procedure for the south¬ 
ern region. This resulted in a solution where the [Neii] 12.8/L(m 
absolute flux was better reproduced (within a factor of 3). This 
model uses an older burst (5Myr), higher density (500 cm“^), 
and the same inner radius (20pc). The corresponding ionization 
parameter is log(t/) = -2.4, which is 0.1 dex lower than pre¬ 
viously found. However, the model underpredicts the neon in¬ 
tensities by a factor of 2 and the [S iv]/[S m] 18.7/um ratio by a 
factor of ~4, while that ratio was matched within 20% without 
the [Nett] constraint. 

5.1.2. PDR lines 

In our default PDR solutions for both regions, [C ii] is system¬ 
atically underpredicted compared to [O i]. The best PDR model, 
found in Sect. I4.2I for the central region single-burst case, in¬ 
cludes a turbulent velocity of 3 km s'. In the southern and cen¬ 
tral region continuous case, similar turbulent velocities also lead 
to better agreement with the [C ii]/[0 1 ] line ratios, but still under¬ 
predict the observed emission in absolute values (AppendixlAli. 

In addition to the stars. X-rays can be a source of heating in 
the PDR and affect the FIR line emission. Point sources and dif- 
fuse X-ray emiss i on ha ve been reported in iHartwell et al.l (l2004h 
and iGhosh et alJ (l2006h . The identified point sources are not co¬ 
incident with the peak of the FIR emission, and we therefore 
ignored them. The diffuse emission is mostly detected in the cen¬ 


tral re gion, with a luminosity of 3 x 10^* ergs * (iHartwell et al.l 
120041) . which is lower than that of the starburst. We have tested 
the effect of this diffuse X-ray component on the PDR lines in 
the central region and found that it increases the predicted inten¬ 
sity of the [O i] lines by ~30% and the [C ii] intensity by less than 
10%. As X-rays are not the main source of heating in the PDR, 
they do not help to produce significantly more [C it] emission. 

We further explored the possibility of [C it] originating from 
a diffuse ionized component. We compared the [C ii] and [N ii] 
122jum intensities and the PACS upper limit on the [N ii] 205/um 
line, which gives [Nii] 122/zm/[Nii] 205/tm >1, to theoret- 
ical predictions assuming pure collisional regime. Following 
iBernard-Salas et aTI (12012h and applying C and N elemental 
abundances observed in NGC 4214, we found that less than 16% 
of the total observed [Cii] emission arises in diffuse ionized 
gas. Herschel SPIRE FTS observations of [N ii] 205//m toward 
the central region also indicate that [Nii] 122jum/[Nii] 205/um 
=2.5 (priv. comm. R. Wu), which gives an ionized gas density of 
~60cm“^ and a contribution of only 8% to the total [Cii] emis¬ 
sion. Therefore, if a low-density ionized component is added to 
our current models to account for the missing [N ii] and [Ne it] 
emission, this component will not contribute significantly to the 
[C It] emission. Our best-fitting models for the central region also 
predict that the H n region contributes less than a few percent to 
the [C It] emission. 

The [C It] emission most likely arises from a neutral phase, 
but its conditions are not well described by our default PDR 
models. With the PDR tests performed in Sect. 14.21 we have ex¬ 
plored the effect of density on the PDR emission lines, but these 
lines are also sensitive to the radiation field strength. To reduce 
the radiation field intensity Go (not the hardness), we placed the 
PDR farther away by stopping the model at the H'^/Hi phase 
transition and resuming the calculation at a larger distance in the 
HI phase (note that this breaks the pressure equilibrium). The di¬ 
rect effect is to dilute the UV f i eld be fore it reaches the PDR, a s 
proposed bv [Israel & Malon^(l201 ih and lCormier et al.l(l20L5l) . 
which is equivalent to increasing the porosity of the medium. 
For the central single-burst central continuous and the southern 
single-burst cases, we increased the PDR distance by a factor of 
2, 3, and 5 (rin=170,186, and 115 pc), respectively. This way, Go 
decreases to ~ 120 in all three cases, boosting the [C ii] emission, 
and the predicted [Cii]/[Oi] ratios match the observed ratios 
within 40% (Fig. |7]l. In absolute values, the [Cii] and [Oi] in¬ 
tensities are 2 to 3 times too high, however. This can be compen- 
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sated by a PDR covering factor lower than unity (see Sect. 15.1.31 
below). 

We conclude that a moderate density, moderate Go neutral 
medium (compared to our dense, high Go default PDR model) 
with a low turbulent velocity and a covering factor lower than 
unity is the most plausible origin for the observed [C ii] emission. 


5.1.3. Model luminosity 

iHermelo et al.l (l2013h have fit the dust SED of the whole galaxy. 
In particular, their models require less UV luminosity than that 
observed to match the IR emission. They discussed possible ex¬ 
planations for this disagreement, proposing an escape of unat¬ 
tenuated UV radiation along with a particular geometry and dust 
properties in the galaxy. Such an argument of a UV escape frac¬ 
tion could also apply in our case. As seen in Fig.|3 the observed 
photometry data points in the FIR, which correspond to emis¬ 
sion originating from the PDR, are lower than the modeled SED 
for the central region. This disagreement can indicate a different 
covering factor for the PDR. In our model, the PDR fully covers 
the sphere around the source (i.e., covering factor of unity). For 
the best-fitting models to better match the photometry, it could 
be that the PDR component is more porous, allowing radiation 
to escape the cloud. To illustrate this, we plot the resulting SEDs 
for the two models of the central region considering a PDR cov¬ 
ering factor of 0.5 (dashed curves in Fig.|6l). 

Part of the discrepancies in our results that we have discussed 
in this section ISTI originate in modeling each complex as a sin¬ 
gle cloud, which is imposed by the lack of spatial resolution 
in the observations. Clearly, future improvements are expected 
from observations with better spatial resolution. 


5.2. Central region: bursty or continuous star formation? 

The star formation history of NGC 4214 over the last Gyr is com¬ 
plex. It shows bursts lasting for shorter or longer periods and 
a continuous ‘back ground’, which takes place throughout thi s 
whole time window iMcOuinn et al.l20T(itlWilliams et al .1201 Ih . 
In that sense, we could say that it is a rather hybrid star formation 
pattern. In the central region, we have investigated cases of both 
a single-burst and a continuous star formation mode. However, 
when we use a single model to reproduce the observed line emis¬ 
sion, we simplify the problem, since we do not take into account 
both modes. Hence arises the question of which of the two ap¬ 
proaches is the most adequate to model the MIR-FIR line emis¬ 
sion. We have shown that both modes can satisfyingly reproduce 
the observed mid- and far-infrared line emission. By comparing 
values of the best-fitting models (Table ID, we found that 
the single-burst case seems to globally perform better when con¬ 
sidering the lines used in the optimization method. When includ¬ 
ing [Nil] and [Neii], the two modes perform similarly, although 
the highx^ values are driven by the poor fit to the [Ne n] 12.8/im 
line. For the PDR, the continuous scenario gives a lower , but 
the default PDR solutions are not optimum and can be fine-tuned 
for both modes (by lowering the density and Go, see Sects. 14.21 
and l5.1b . We conclude that both are limiting, simplifying cases of 
modeling the ISM in NGC 4214, with the continuous star forma¬ 
tion model being marginally more accurate inside the PDR and 
the single-burst in the H ii region (without further refinement). 


Table 7. Comparison of star formation properties. 


Quantity 

Region I 

Region II 

Tcii/Ttir 

5.4 X 10-3 

3.6 X 10-3 

Tcii/Tco(i- 0 ) 

6.7 X 10"^ 

2.5 X 10^^ 

Mgas=m+HI [Mq] 

2.05 X 10® 

2.68 X 10® 

Mh,/Mhi 

0.35 

0.32 

SFR [Moyr-'l 

2.2 X 10-3 

1.9 X 10-3 

SFE [Gyr-'l 

10.7 

7.1 


5.3. Comparison between the two star-forming regions 

How do the ISM conditions that we have characterized in the 
two star-forming regions relate to their star formation proper¬ 
ties? We have found that the modeled cluster (radiation source) 
in the southern region contains younger stellar populations with 
a harder radiation field than that in the central region, in agree¬ 
ment with the results of llJbeda et al.l (12007 ah . for instance. The 
hydrogen density is also higher in the southern star-forming re¬ 
gion, but the metallicities of the regions are very similar. The 
southern region is observed to be at a younger, more compact 
stage than the central region. The central region is more evolved 
and had time to expand, as observed by the p resence of shells tha t 
may have swept away the dense material ([Walter et al.l 1200 ih . 
and is thus consistent with a more diffuse ISM. 

We calculated the star formation rate surface densities for the 
two regions combinin g the GAFEX FU V map and the Spitzer 
24jum map, as done in Ferov et aH (l2008h . with 

5:sfr[Mo yr-' kpc-^] = 3.2 • lO'^ • I 24 [MJy sr-‘] 

-H 8.1 ■ 10^2- Ifuv [MJysr-'l, (3) 


where Ssfr is the SFR surface density and I 24 and Ifuv the 24/im 
and FUV intensities. We also measured the atomic and molec¬ 
ular hydrogen co ntent of the two re gions using the 21 cm map 
from THINGSEI ([Walter et al.l 1200^ and th e CO(l-O), CO( 2 - 
1) transition maps from [Walter et al.l (1200111 and HERACFEl^ 
(iFerov et al.[[2009[l . respectively. We used a conversion factor of 
ffco=4.38 [Mopc^'(Kkms"*)"'] from CO(l-O) luminosity to 
H 2 mass. If we were to use a different co nversion factor due t o 
the low metallicity of these regions (e.g., ISchruba et al.[[20l'^ . 
this would not affect the relative comparison of the regions (see 
Table |7]i. The central region (I) has a total hydrogen content of 
Mgas,i=MHi-(-MH 2 = 2 .05x10® Mq and a molecular (mass) fraction 
of /inoi,i=MH,/MHi=0.35. Integrating Eq.|3]in the region, we find 
SFRi=2.2xlO“^MQyr '. The southern region (II) has a higher 
total hydrogen content of Mgas,ii=2.68x10® Mq, an H 2 fraction 
/moi,ii=0.32, and SFRn=1.9xlO“^Moyr '. Therefore the south¬ 
ern region has relatively more gas compared to its SFR than the 
central region. In terms of efficiency, SFR/Mgas, it is about 50% 
lower in the southern region. This could reflect a slightly more 
efficient, cluster-like star formation episode in the central region 
or simply encode a different evolutionary state as the SFR and 
gas masses are a strong f unction of time on sc ales of individual 
star-forming regions ('e.g.. [Schruba et al.l2010[) . The southern re¬ 
gion being younger, it may still be in the process of forming 
stars. 

The main differences in ISM conditions that we extracted 
from our modeling relate to the Hii region properties. The emis¬ 
sion lines are a factor of two lower in luminosity in the south¬ 
ern region, except [Neiii] and [Siv], which are proportionally 


^ http: //WWW. mpia- hd. mpg. de/THINGS/Data. html 
^ http://WWW.cv.nrao.edu/~aleroy/heracles_data/ 
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Central region - burst Central region - continuous South region - burst 



Observed 1^^ Predicted (default) Predicted, D(PDR) increased 


Fig. 7. Effect of changing the distance to the modeled cloud for the PDR calculation. The PDR distance is increased hy a factor of 2, 3, and 5 for 
the central single-burst (left panel), central continuous (middle panel), and southern single-burst cases (right panel), respectively. 


higher in the southern region. By contrast, the PDR conditions 
in the two regions are similar (hh - 10"*cm“^ and Go - 150). 
Our modeling reflects conditions resulting from the recent star¬ 
forming event and has little predictive power regarding a differ¬ 
ent, future star-forming event. In particular, at the linear scale 
that we probe (~175pc), the PDR conditions are averaged over 
multiple star-forming clouds and not representative of the un¬ 
derlying, possibly different, substructure in individual molecular 
clouds. However, there is more PDR emission relative to CO 
in the ce ntral region, that is, high [C iil/C O and [Oil/CO ratio s 
(see also lCormier et alJl201(ll) . As found by IWalter et al.l (1200 ll) . 
CO emission is centrally concentrated in the south and more dif¬ 
fuse in the center. The concentration of molecular gas may be 
nourishing the current star formation episode in the south or is 
being observed at a pre-disruption stage with the same fate as 
the central region. The increased porosity, which evidently is an 
intrinsic property of the low-metallicity ISM, is seen in both the 
central and southern star-forming regions. The main evolution 
within the dense medium is seen in the covering factor of the 
PDR, which is found to be lower in the central, more evolved 
region than in the southern region, and in CO, which probably 
suffers more from photodissociation with time and its emission 
is seen farther away from the cluster center, but this cannot be 
modeled with our static approach. Observing an intermediate 
PDR tracer at the C^/CO transition, such as C i, would help to 
test this evolution scenario. 


6. Conclusion 

We have investigated the physical conditions characterizing the 
ISM of the dwarf irregular galaxy NGC42I4 by modeling 
Spitzer and Herschel observations of MIR and FIR fine-structure 
cooling lines. We used the spectral synthesis code Cloudy to 
self-consistently model the H ii region and PDR properties of the 
two main star-forming regions in NGC 4214. We summarize our 
results as follows: 

- The ionized gas in the southern region is found to be 2.5 
times denser than in the central region (440 cm“^ versus 
170cm“^) and typified by a harder radiation field. Our best¬ 
fitting models of the H ii region+PDR reproduce most ionic 
and neutral atomic lines, namely the [Out] 88//m, [S iii] 18.7 
and 33.5pm, [Siv] 10.5pm, [Neiii] 15.6pm, and [Oi] 63pm 
lines, within a factor of ~2. 


- The observed [Nii] 122pm and [Neii] 12.8pm lines are the 
most discrepant with our model solutions for the H ii region. 
A single model component seems too simplistic to account 
for all observed lines simultaneously. Given the complexity 
of these star-forming regions, a multi-component modeling 
would be more appropriate. In particular, a lower excitation 
ionized gas component may be required to match the [N ii] 
and [Ne ii] emission in both regions. 

- Our Hit region models and the established observational 
[C ii]/[N It] line ratio used as a proxy for the fraction of [C ii] 
arising in the ionized gas both indicate that the [C ii] emission 
is mostly associated with the PDR gas, with only a ~10% 
contribution from the ionized gas. 

- Constant pressure models where thermal pressure dominates 
the pressure equilibrium perform rather poorly for the PDR 
lines mostly because of the high densities and high Go val¬ 
ues reached in the PDR. Including additional pressure terms, 
such as weak turbulent or magnetic pressure, or placing the 
PDR cloud farther away and reducing its covering factor, 
leads to a much improved reproduction of the observed line 
intensities. 

- Star formation histories have an effect on the predicted MIR- 
FIR line emission. We have explored the two simplifying 
cases of a bursty and a continuous star formation scenario. In 
the central region, we found that the bursty scenario works 
marginally better for the H ii region and the continuous sce¬ 
nario for the PDR, although both modes can reproduce the 
observations after refining the PDR conditions. 

- The H II region modeling from IR emission is consistent with 
the evolutionary stages of the regions found in previous stud¬ 
ies: the southern region is younger and more compact, while 
the central region is more evolved and diffuse. On the linear 
scale of our study (~175pc), the PDR conditions of indi¬ 
vidual star-forming clouds are averaged out and do not echo 
the observed differences between the two regions (stellar 
ages, Hii conditions, etc.). The increased porosity of the star¬ 
forming regions appears as an intrinsic characteristic of the 
low-metallicity ISM, while the covering factor of the PDR, 
which is reduced in the central region, stands out as the main 
evolution tracer. 
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Appendix A: Comparison of PDR modei predictions 
and observations 

We show in this appendix the effect of turbulence and magnetic 
fields on the predicted PDR line intensities for the central re¬ 
gion continuous case and for the southern region single-burst 
case. Tests with turbulence and magnetic fields are presented in 
Sect. l4.2l onlv for the central region single-burst case. 


Central region - continuous 



South region - burst 



Fig. A.l. Same as Fig. [5] for the central region continuous case (top 
panel) and for the southern region single-burst (bottom panel). See 
Sect.l4.2lfor details. 
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